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(Dated: February 5, 2008) 

Spontaneous switching events in most characterized genetic switches are rare, resulting in ex- 
tremely stable epigenetic properties. We show how simple arguments lead to theories of the rate 
of such events much like the absolute rate theory of chemical reactions corrected by a transmis- 
sion factor. Both the probability of the rare cellular states that allow epigenetic escape, and the 
transmission factor, depend on the rates of DNA binding and unbinding events and on the rates of 
protein synthesis and degradation. Different mechanisms of escape from the stable attractors occur 
in the nonadiabatic, weakly adiabatic and strictly adiabatic regimes, characterized by the relative 
values of those input rates. 



Information may be passed from one cellular gener- 
ation to another not just in the form of the DNA se- 
quence, but also in the long lived expression patterns of 
genes. The epigenetic state of the cell i.e. which genes 
are expressed at a given time, is determined in part by 
binding and unbinding of transcription factor proteins 
to the DNA. The genes with their partner proteins form 
complex dynamical systems known as genetic networks, 
which can have many steady states i.e. an attractor 
landscapei*^. The attractors are more stable than the 
individual molecular protein-DNA adducts, because the 
protcomic atmosphere of gene products renews the DNAs 
binding state, ultimately creating auto-catalytically its 
own proteomic atmospherai^. The attractors of such 
a genetic network may be associated with distinct cell 
typea^. Experimental evidence for this view has re- 
cently been presentedSii. The growing experimental in- 
terest in this problem^, as well as a number of theoretical 
puzzles involving the stability of the attractors^'^, call for 
a flexible and intuitive theory of the lifetime of such ge- 
netic network attractors. Some progress has already been 
made towards the goal2ii2iii*iS, but existing formalisms 
are cumbersome, certainly when compared with the the- 
ory of activated events in molecular systems based ulti- 
mately on transition state ideas Our goal here is to 
present a simple treatment of the noise induced transi- 
tions between two attractors on a landscape that is par- 
allel to the treatment of simple molecular rate processes, 
which starts with Wigner's absolute rate theor}ii^. In 
chemical kinetics, the ratio of escape is proportional to 
the probability of rare configurations equally likely to 
become reactant or product. These rare configurations 
represent a stochastic separatrix of the motion. 
While thermal atomic motions cause the escape from 
energy minima in molecular physics, the noise in ge- 
netic networks comes from the probabilistic nature of the 
chemical reactions, since only a small number of proteins 
and individual copies of the target DNA are involved. 
Unlike molecules, genetic systems being far from equilib- 
rium, cannot strictly be described by thermodynamic free 



energy functions. The stochastic separatrix for molecu- 
lar activated events is a dividing surface passing through 
saddle regions of the free energy. We argue that, even 
in the absence of a free energy function, the notion of 
a stochastic separatrix between basins of attraction re- 
mains a good approximationi"'^^ and allows a treatment 
of stochastic switching along the lines of a transition state 
theory with dynamic corrections involving the rates of el- 
ementary processesl2ilSi. 

The dynamics of gene networks involves two very dif- 
ferent processes whose rates must be compared- protein 
synthesis and DNA binding. The complexity and en- 
ergy consuming nature of protein synthesis, in prokary- 
otic cells, generally causes changes in protein number to 
take longer than the diffusion controlled binding times of 
transcription factors, even at their low concentrations. 
For this reason, it has been argued that one can de- 
scribe the binding of the transcription factors to each 
DNA binding site as an instantaneously equilibrated pro- 
cess, when considering protein production. For steady 
states this approximation appears to be reasonably accu- 
rate. It has, however, also been notediiiSii2iS£ , that the 
DNA state fluctuations may influence the protein num- 
ber state fluctuations. Here we show that the impact 
of DNA state fluctuations on the escape process, is con- 
siderable in a rather wide parameter regime for which 
the steady states are not much influenced by the DNA 
state fluctuations (Figure . For biologically relevant 
parameters, the DNA occupancy fluctuations may sig- 
niflcantly increase the spontaneous switching rate from a 
given attractor. In what we call the nonadiabatic limit, 
where DNA state fluctuations dominate the protein num- 
ber fluctuations, individual binding and unbinding events 
of the transcription factors are directly responsible for the 
transition. For much of the adiabatic regime, although 
the influence of DNA fluctuations on the steady state pro- 
tein levels is negligible, these fluctuations still modify the 
lifetime of a state - we will call these transitions "weakly 
adiabatic" . DNA occupancy fluctuations can only be ne- 
glected at very high values of the rate ratios, in what we 
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FIG. 1: The sum of the escape rates k — kon + koff as a 

hg'^ 

function of the adiabiaticity parameter k = for a self- 
activating switch with g-\ = 100, (/| = 8, fc = l,n]y = 53.4. 
Comparison of the exact discrete n numerical calculation 
based on the mean free passage time (black solid line), with 
approximate methods: in the nonadiabatic limit (small k) 
(gray dashed line, Eqs Q and |3Jl, in the weakly adiabatic 
regime (black dashed line, Eqs |21 and (HJ and mixed crossover 
regime (gray solid line). The adiabatic results tend asymptot- 
ically to the strictly adiabatic limit (large K-flat escape rate) 
(light gray dashed line, Eqs|S|and|HJ. In the strictly adiabatic 
limit the binding of transcription factors to the DNA binding 
site is equilibrated. In the nonadiabatic and weakly adiabatic 
limits the escape rates show a dependence on the adiabatic- 
ity parameter- the process in influenced by the DNA binding 
state fluctuations. 



call the strongly adiabatic limit. 

As Sasai and Wolynesi have pointed out, the stochastic 
theory of a simple genetic switch, can be considered anal- 
ogous to the physicists' Kondo problem or the chemists' 
electron transfer process, where DNA occupancy plays 
the role of a spin or electronic state variablei^. In a sim- 
ple and intuitive way, here we exploit these analogies to 
compute the lifetime of a genetic switch, using the idea of 
a landscape with a stochastic separatrisiS, much like in 
the earlier proposed threshold modeli^. Our treatment is 
quite analogous to that used for characterizing adiabatic 
vs nonadiabatic regimes of quantum ratesi^^i^. The ap- 
proach, we present is easily generalizable to a many gene 
system. In the general case the present approximation 
yields the lifetime of a given state of the switch, which is 
governed purely by a few local properties of the landscape 
and does not require computing complicated trajecto- 
ries. Global properties, sensed by the most probable es- 
cape paths, generally enter rates for far-from-equilibrium 
systems^. The present approach must, therefore, be ad- 
mitted to be approximate. The simplicity hopefully will 
make up for some inaccuracy. 

Several treatments of the mean first passage time be- 
tween epigenetic states have already appeared. Most of 



these studies assume the DNA state equilibrates on a 
much faster time scale than the protein numbes^SiSi. We 
refer to this as the adiabatic limit. In this limit the pro- 
tein number states may then be treated as a continous 
variable giving an expression for the mean first passage 
time a la the Smoluchowski theory of diffusive rates as 
sketched by Bialek-'^^. A more rigorous approach finds 
the rate by constructing the most probable escape pathS 
or by calculating the distribution of pathsiSi^. These 
methods are powerful, but they are hard to visualize, 
especially for more complex switching systems. While 
the usually invoked adiabatic limit seems to be appropri- 
ate for simple switches in prokaryotes, it is not an obvi- 
ously correct approximation for switches that have more 
complex operators, in which multiple protein elements 
must combinatorically assemble at a given site, slowing 
the binding^^. Nonadiabatic effects should also play a 
significant role in eukaryotic systems where chromosome 
restructuring, which may be quite slow, dominates the 
epigenetic transition. Artificially engineered switches^^ 
may be constructed with parameters spanning the entire 
phase diagram. 
The Simplest Switch 

For illustration we will present our ideas using the sim- 
plest example of a system in which we can consider the 
escape from one minimum to another - a bistable self- 
activating switchSfli. We emphasize the approach is more 
generally applicable. The self-activating switch consists 
of a single gene, which may be found in one of two states: 
on or off. In the off state proteins are produced at a basal 
level, but in the on state proteins are produced at an en- 
hanced level, leading to a number, n, of proteins in the 
cell at any moment. The proteins act as activators by 
binding to the same operator site as the gene governing 
their production. We assume they bind as dimers with 
a rate h(n) — hn{n — l)/2. The unbinding of the tran- 
scription factors is described by a rate /. We neglect 
time delays due to mRNA synthesis etc. (which admit- 
tedly may play a key role), so that protein population 
dynamics is governed by a birth death process. Protein 
degradation occurs with a rate k, production with the ac- 
tivated rate gi; in the on state and the basal rate gi in the 
off state. The system is characterized by a two state joint 
probability distribution P{n), describing the probability 
of having n proteins in the system and the DNA binding 
site being in the bound (on-t) or unbound (off-|) state. 
A recent combined experimental and theoretical study2& 
has brought attention to the bistability of a switch in 
previously unexplored limits, when the degree of operon 
repression is small. Our discussion will turn also to the 
nonadiabatic limit. Here the equilibration of the DNA 
and changes in the protein number occur on comparable 
time scales. 

To compute escape rates from the steady state attractors 
one must determine the stochastic separatrix^^.. In the 
adiabatic limit, the position nj^ of the minimum of the to- 
tal probability distribution P(n) = P^{n)+Pi{n) is given 
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by the condition of zero mean protein flow dn/ dt\\ ^ = 

(M + h{n)g^)l{f + h{n)) - kn\^^^^^ = 0. For a bistable 
switch, this equation is satisfied by three values of n; - 
one solution gives the separatrix, the other two the po- 
sitions of the high and low protein number stable steady 
state attractors, and n\. In the nonadiabatic limit 
the stochastic separatrix refers both to the DNA and pro- 
tein number state. This results in a different value of the 
critical separatrix numbers njy in the nonadiabatic, and 
n\ in the adiabatic limits. The direction of flow changes 
when the DNA state changes. Therefore the position 
of njy corresponds to that number of proteins needed 
for the system to have comparable probability to be in 
the on or the off state. For simplicity we can approxi- 
mate in the large n limit h(n) = h/2n(n — 1) « hii? /2 
and determine the position of the nonadiabatic separatrix 
by means of mass action, using the chemical equilibrium 
constant K'"!: nj^ = VVk^, where K'"i = 2f/{hV^), 
where V is the cell volume. The steady state attractors 
in the nonadiabatic limit are determined by the birth- 
death processes in the particular DNA states: = 9i/k 
in the off state and n'^ — /k in the on state. To func- 
tion as a switch must be less than njy and must be 

greater than njy. We can rewrite the adiabatic separa- 
trix positions in terms of the volume scaled equilibrium 
constant K^'^V^, which scales with nj^, as n\ — K'^'^ j-n) 
and n-^ < < < njy < < n^. 
Nonadiabatic Rate Theory 

Here we compute the rate of escape of the system from 
the low protein number attractor to the high protein 
number attractor ikon) and vice versa (koff)- Since in 
the nonadiabatic limit the low protein number attractor 
corresponds to the off DNA occupancy state and the high 
protein number state corresponds to the on DNA occu- 
pancy state, the transition from the low protein number 
state to the high protein number state is requires the 
binding of an activator. Without the possibility of bind- 
ing and unbinding, the dynamics in each attractor would 
be described by stochastic destruction and production of 
proteins alone, resulting in fluctuations of the mean pro- 
tein number around each steady state. Consider a sys- 
tem maintained in the off DNA binding state and that 
now has proteins. The initial probability of being in 
the off DNA state, with precisely proteins present is 

-Po/ /("''^) = -Pi('^'^)/(^T('*'^)~'"^i("''^))- "^"^ may be gener- 
ally assumed to be close to the mean number of proteins 
in the off state (n^ = 9l/k)- If a binding event now 
occurs at time t — 0, the gene spontaneously flips into 
the on state and proteins are now produced at an en- 
hanced rate. The protein number increases towards the 
mean number in the high protein state (n^ = gi/k). If 
the activator does not unbind before the number of pro- 
teins becomes characteristic of the on state attractor a 
successful switching event will have taken place and the 
protein number will now fluctuate around the on steady 
state value. However, since, we are in the nonadiabatic 



limit, the timescales to reach the steady state for both 
the DNA binding state and protein synthesis and degra- 
dation are assumed comparable, so an activator may in 
fact unbind before reaching the separatrix at n]y. If an 
activator does unbind during that time, the gene returns 
to an off state, albeit with a slightly higher number of 
proteins than initially. Another binding event will re- 
peat the above scenario, until the protein number safely 
crosses the separatrix at rijy and the steady state corre- 
sponding to an activated gene is reached (Figure [5] a) . 
The average time needed to cross the barrier from an 
initial point n-^ , which is also the time allowed for a un- 
binding event to occur, is the mean time to reach nJy for 
the enhanced production rate. The initial rate of bind- 
ing an activator /i(n^) = h/2n^{n^ — 1) must be modified 
to account for the possibility of unbinding again before 
the system crosses the separatrix. Summing of these at- 
tempted crossings, results in an expression for the rate of 
escape from the off state minimum (n^) to the on state 
(n > njy) in the nonadiabatic regime given by: 

fco„(n^)=Po//(n^)/i(n^)e'J'("" (1) 

The exponential term gives the successful fraction of 
attempts to reach the protein number based separatrix, 
launched from the steady state The total time to 
reach the separatrix is given by i(?ijv) ~ as deter- 

mined by the average flows in the initial DNA state and 
the mean time for an unbinding event to occur is /~^. 
Explicitly, the escape rate from the off state, becomes 
kon{n^) = Poffin^Mn^KigT - kn^)/{gT - kn^^))-f/K 
The power-law term describes the motion on the surface 
with enhanced production after binding of the activator. 
In the nonadiabatic limit, the probability distributions 
for the on and off states are unimodal. Therefore it is 
unlikely for the gene to be in the on state is the number 
of proteins is small, thus poff{n-^) ~ I. If the protein 
number is large and the unbinding rate is comparable to 
the death rate this expression yields: 



where k — hg'^ / (2k^) . In the extreme nonadiabatic limit 
K — > 0, the first attempt may be successful hence the 
result simplifies to fco„(n^) ~ h{n^). 
A similar calculation can be carried out starting from the 
other steady state. The escape rate from the on state, 
with > njy proteins, is given by the rate of binding of 
an activator at time t = 0, providing the system is in the 
on state Pon{n') — Pt {fT-^ / {P^ {^^ ) + Al"-^))' reduced 
by the probability that an activator rebinds before the 
protein number decreases to numbers characteristic of in 
the off state {n < njy). The time available to rebind is 
calculated using protein production at a basal level. The 
koff rate is therefore: 
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For the off rate the mean free path before a rebinding 
event depends on the mean number of proteins in the 
system n. The argument of the exponential still de- 
scribes the number of rebinding events. In the strongly 
nonadiabatic case, Pon{n^) ~ 1, and for very large mean 
protein numbers the escape rate tends to: 

koff{n^) - /e-3^("' -< ) ^ /e"^— TT^ (4) 

Due to the timescale separation in the nonadiabatic limit 
the system may be approximated as a two state system. 
The ratio of the escape rates, therefore yields the ra- 
tio of the probabilities to be in the individual steady 
states. The equilibrium constant for the "dressed" ge- 
netic states in the nonadiabatic limit K^^ = kofj/kon 
therefore becomes K'~^^ {n\^/n-^)'^exp{—K/2) = 

K^'^V'^ /n^'^exp{~K/2). When k = the proteomic at- 
mosphere has no effect on the relative stability of the 
DNA occupancy, which follows the ordinary mass action 
law. 

The formulae described above provide quite intuitive rep- 
resentations of specific escape mechanisms. These results 
may also be formally obtained via the path integral so- 
lution of the master equation by using the method de- 
scribed by Wang, Onuchic and WolynesSi for kinetic pro- 
tein folding. This result also coincides with the heuristic 
approach of Ninic'^*. 

Adiabatic Rate Theories: Weak and Strong 
Regimes 

In the nonadiabatic limit the switch reaches the separa- 
trix within the time for a few binding events, as schemat- 
ically portrayed in panel a of Figure El In what we 
call the weakly adiabatic regime, the escape process pro- 
ceeds differently. The DNA occupancy responds quickly 
to the changing proteomic atmosphere reaching a local 
steady state before the protein number changes by a 
large amount. The average occupancy then determines 
the average local rate of protein synthesis and degrada- 
tion. A few binding and unbinding events are required in 
the nonadiabatic limit, but in the adiabatic limit those 
events are much too common to allow the direct mecha- 
nism. One is tempted to equate the local diffusion rate 
to that coming from synthesis and degradation. But this 
temptation can only be rigorously indulged at an extraor- 
dinary high binding rate. Instead a random, but cyclic 
process of binding, growth and unbinding churns the pro- 
tein number like a turbulent surf. The cyclic motions of 
eddies in an ocean wave, if interrupted contribute to a 
diffusive transport of flotsam to the shore. In the same 
way, in most of the weak adiabatic regime, protein num- 
bers fluctuate from the mean flow through this "churn- 
ing mechanism" . The protein number, changes slightly 
with each cycle of binding/growth/unbinding and even- 
tually reaches the separatrix point due to the resulting 
diffusive motion. One can show the system acts as if it 
were diffusing along an effective potential, whose gradient 
gives the mean flow expected from the average occupancy 
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FIG. 2: A schematic diagram of the difference in the charac- 
ter of the transitions from the state with a small number of 
proteins to the state with a large mean steady state number 
of proteins in the nonadiabatic (A) where the escape rate is 
given by EqsQand|3 adiabatic (B), where the escape rate is 
given by Eqs|^and|Hl and extremely adiabatic (C) regimes. 
The dark gray line marks the effective potential for protein 
number change. The horizontal arrows signify binding and 
unbinding events. 



V{n) — geffin) — kn (panel b in Figure|21). The diffusion 
rate in this outwardly adiabatic regime though depends 
on the nonadiabatic events. Only at very high adiabatic- 
ity is diffusion ascribable to birth-death alone. 

It is helpful to understand the "eddy-induced" dif- 
fusion in an intuitive way. The effective production 
rate geff{n) = {fgi + h{n)g^)/{f + h{n)) is the pro- 
duction rate averaged over the binding and unbinding 
states, if they were in equilibrium. The diffusion ex- 
pected solely from the birth-death processes would just 
be DBD{n) = geff(n) + kn. This fluctuation mecha- 
nism is augmented by diffusion in the orthogonal two 
state "binding-space", that is the eddy motion. The 
difference in the mean distance in protein number that 
would be travelled in the two DNA states during a typ- 
ical eddy cycle will be An = (g^ — gi)/{f + h{n)). 
It is the typical difference in protein number expected 
after a full cycle of an eddy has been traversed. It 
is given by the difference in velocity in protein num- 
ber space in a given binding state, Au = |5| — 
times the mean time before the binding state changes 
M = {f + h{n))-^, such that An = AwAi. The mean 
free time, or the eddy mixing time, is given by the sum of 
the characteristic times for binding and unbinding, both 
of which must occur to return to the original binding 
state, T = /^^ -I- {h{n))~^ . The rate which describes the 
eddy cycling thus becomes — fh{n)/{f + h{n)). The 
diffusion coefficient D = An^/r is the square of the mean 
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change in protein number divided by the characteristic 
time spent within a given eddy. The latter depends on 
both binding and unbinding events. One thus finds in 
DMndzngin) = fh{n){g^ - gifl{f + h{n)f . 
The mean number of proteins of a given type produced 
in the active state is of the order of g^/k ^ 10^. The 
degradation rate of proteins gives hfetimes of the or- 
der of a bacterial generation k ^ 10^^ s^^. Dissocia- 
tion rates from the DNA vary from / 1 — 
and typical equilibrium constants may be taken to be 
K'^iy'^ ~ 10^ — 10^, which results in association con- 
stants h/2 = f/{K^W^) ~ 10-2 - (based on A 
phage data as assembled inSS and references therein). We 
therefore see that typical adiabaticity parameters scan a 
wide range: k = hg'^/{2k^) ^ 10° — 10^. The diffusion co- 
efficient from churning, which depends on the DNA occu- 
pancy dynamics, typically influences the escape rate over 
four orders of magnitude of the adiabaticity parameter 
K G (lO'^ — 10^), nearly covering the biologically relevant 
regime. For escape processes the DNA binding dynam- 
ics cannot be neglected until the adiabaticity parameter 
becomes extremely large ultimately yielding the strongly 
adiabatic regime. As shown in Figure |21 the eddies due 
to the influence of the DNA binding state become smaller 
with faster binding, until the motion becomes dominated 
by simple birth-death diffusion along the effective poten- 
tial, giving the steady state probabilities, averaged over 
the DNA binding states (panel c of Figure ^ . 
In the adiabatic limit, the escape rate is governed largely 
by the fraction of systems at the separatrix Nj^ com- 
pared to the fraction residing near the original attrac- 
tor Nattr- N\/Naur = P(n^ )/p^< (n^) , where p"<(n^) = 
Sn<nt -^C^^) ^"^^ is the steady state probability 

density for a state with n proteins. Clearly p<('^^) = 
P{nin)5nin, where Suin is the width of the attractor. 
It is important to understand the spatial variation of 
P(n), described by the "potentials" in Fig. [3 The spa- 
tial variation depends on the balance of the local mean 
flow against the flow due to diffusion. We can under- 
stand this balance by considering the motion pictured in 
Figure |21 h. The mean local velocity by which the pro- 
tein number changes is u = ffe// — kn. In addition to 
this drift the protein number changes by diffusive mo- 
tion from places of low to high probability, with a ve- 
locity of Vdif fusion — '2D^(n)/lc, where Ic is a charac- 
teristic "lengthscale" over which the steady state prob- 
ability changes by roughly one e-fold. D^{n) refers to 
the diffusion coefficient, which governs the motion in 
a particular regime. It is equal to Dunding{n) in the 
weakly adiabatic regime, DBoin) in the strictly adia- 
batic regime and is roughly DsDin) + Di,inding{n) in the 
small crossover region in between. To traverse this scale 
the local velocity has to be at least as large as the ve- 
locity of the diffusive motion v > Vdif fusion- The equal- 
ity V = Vdif fusion sets a characteristic length scale of 
the problem Ic = 2L)'(n)/|w(n)|, over which locally the 
probability in a steady state should change by a factor 



of e. This relation is valid both in the adiabatic and 
nonadiabatic regimes. The quantity Ic is analogous to 
the "scale height" in the equilibrium barometric prob- 
lem. How many of these characteristic steps of length Ic 
are needed in order for the system to reach from its 
steady state value? Bearing in mind that the length of 
each step, depends on n, we must concatenate these steps 
to give the probability to be at the separatrix relative to 
being near the initial state. The probability exponen- 
tially depends on the number of scale heights of varying 
length Ic needed to reach the improbable separatrix start- 
ing from the most probable situation at the basin center, 

exp[-2/;/ 

"■A 

To find the rate, we finally need the width 5nin. The 
size of the attractor Sriin is analogous to Ic at the bot- 
tom of the basin, but quadratic order effects must be in- 
cluded. To compare the velocities of the motion near the 
basin center due to drift and diffusion, the drift velocity 
must be computed as the "drift frequency" in the initial 
state w(n„) = {dv{n) / dn)\n,„ = fhnin{g^ - 9l)/if + 
h{nin)Y — k times the distance from the stationary point. 
Comparing drift and diffusion velocities in the same re- 
gion uj{nin)5nin = D''(nin)/Snin, gives the size of the 
attractor Sriin = \/\D^ {n-in ) / w (ni„ ) | . The exponential 
term counts the paths from all possible position within 
the attractor. We must therefore divide the by the width 
of the attractor. 

To determine the epigenetic escape rate we need also 
the transmission factor. In the adiabatic limit, reach- 
ing the separatrix does not yet guarantee a successful es- 
cape. Once the protein number reaches the vicinity of the 
stochastic separatrix the system may directly cross the 
separatrix, or recross it many times before committing to 
the new attractor. The number of escapes per unit time 
rate is thus proportional to the velocity with which the 
system moves over the separatrix, divided by the number 
of attempts before it successfully commits to the new at- 
tractor k — 5v/NP{nin n peak). The velocity 
around the peak is determined by a mean free path for 
number fluctuations Imfp and a mean free time r relevant 
to that region, 5v — Imfp/T. Only in the crossover region 
is it necessary to take all processes into account on equal 
footing when evaluating the mean free path Imfp and the 
associated mean free time r. In the weakly and strongly 
adiabatic limits the results simplify. In the weakly adi- 
abatic region, the mean free path is dominated by the 
DNA churning cycles and is given by the typical eddy 
size l„,fpK. {g^-g^)/[f + h(n)) and r = /"^ -h (/i(n))-^ 
In the strictly adiabatic limit, the motion is determined 
by the birth and death of proteins. Effectively, the pro- 
tein number changes by Imfp equal to one protein in the 
mean free time r = {geff{n) + fc)"^. Once the mean 
free path has been determined, the number of crossings is 
then the number of steps of the size of the mean free path 
needed to cross the transition state region Itst- Like the 
basin size, the size of the transition state region is Itst = 

\J D^{n\)/Lu{n\). The escape rate from the left attrac- 
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tor, n^, is kon 



llfp/{lTSTr){5n-\,) 



-n 



dnl 



where I'^fp/''' ~ ^^C^-a) ^^"^ * indicates BD, binding 
and mixed in the appropriate regimes. This gives the 
rate of the escape from the low protein number state in 
the weakly adiabatic regime: 



konin-A)^T^D'{n^^)., 

2vr y D^{n\)D^{n\ 



e ^ 



(5) 

where is the number of proteins corresponding to the 
minimum of the total steady state probability distribu- 
tion. In the adiabatic regime the separatrix is given as 
the fixed point of the average flow: geffin]^) ~ kn\. 
In the strictly adiabatic limit, the eddy motion may 
be neglected. So is determined solely by the 

equilibrated diffusion in protein number space 1'^^°° k, 
ideff + kn)/{geff — kn) . All the components in Eq[Sl 
can be obtained using quadrature, in this case, yielding 
a complex expansion. A more simplified result, explicit 
in terms of chemical rate constants, follows if we linearize 
ZjT^ in the region, which contributes most to the result of 
the integral. In this situation equation [S] becomes: 
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FIG. 3: A phase diagram as a function of the activated pro- 
duction rate g-\ and the unbinding rate / for constant K'^'^V^ , 
showing the areas of parameter space where a given escape 
mechanism dominates based on the ratio of the size of the 
transition state region Itst to the mean free path Imfp- If 
the number of crossings of the separatrix is large i^isx- > 
the transition is adiabatic. If the system commits to a new 
attractor after one crossing '^-^^ < 1 the transition is nona- 
diabatic. 



(6) 



where nmin number of proteins for which l'^ has 
the largest value. The largest value of cor- 
responds to the the smallest characteristic length 
scale in the region of integration. The value 
of scales as ~ V^K'"i/2. The 

pre-exponential factor has the form fi{K^'^) ~ 

kV/{2TT) ^J K<'i/{aon'^\n'^^ - {K'^^V'^Y - 2i^«9F2(nT)2), 
where oq = gilg\- The escape rate decreases with the 
equilibrium constant and system size. Using the depen- 
dence of the minimum of the integrand as a function of 
the equilibrium constant K^'^V'^ ^ one finds the escape 
rate scales as e^"^"^ (if v^-Saon^^ )^^^^ where a\ is a 
numerical factor of the order of 1/2. The rate of escaping 
from the off state attractor exponentially decreases with 
increasing of the equilibrium constant. 
How the escape rate depends on the molecular param- 
eters, can be seen by assuming, for simplicity, a highly 
cooperative variation of the equilibrium DNA occupancy 
with protein concentration. In this case the effective 
production rate can be approximated by the production 
rate in the off state attractor, g'^^^in) k, gi^. Now, 
the protein dynamics will be determined by the rates 
characteristic of the attractors, until the system reaches 
the separatrix. This approximation is like the threshold 
picture of Metzler and Wolynesi^. In this approximation 
one finds: 



A' 



(7) 



When the cell is sufficiently small the separatrix merges 
with both attractors. In such a regime, this simple for- 
mula correctly predicts the functional dependence of the 
escape rate on the equilibrium constant and the pro- 
tein production rates. When the separatrix begins to 
merge the attractor, the exponential term approaches 
unity. Thus stability is compromised. When the attrac- 
tors merge with the separatrix the pre-exponential factor 
becomes important for quantitative analysi a^i^^ . 
In the K dependent weakly adiabatic region, the prob- 
ability distributions look qualitatively similar to those 
in the strictly adiabatic limit: the extrema do not 
change as n increases. In the escape rate calcula- 
tion, however, one compares the ratios of the prob- 
abilities near the minimum and the saddle regions. 
This ratio is significantly different in the weak and 
strong adiabatic regimes and strongly affects the spon- 
taneous switching rates, as seen in Figure ^ In 
the weak adiabatic regime one finds the escape rates 
depend exponentially on the adiabaticity parameter 
K. The escape rate therefore is approximately dom- 
inated by kV/{2'K)K'"in^^/aay/n^^ - 2K<'iV^ / {n"^^ + 
K^W^)^ expi-f /kK^iV^ao/n^n\ - n^^) / {n^^n\)) . In 
the weakly adiabatic regime, the effective growth rate 
can be well approximated as that with a fixed DNA oc- 
cupancy, as in the Metzler- Wolynes threshhold modeU^. 

The transition can be treated from the high protein 
number state to the low protein number state much as 
above in the adiabatic limit. The rate of escape from high 
protein number to the low protein number depends on 
the relative probability that the system is to the right of 
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the separatrix, characterized by a mean protein number 
compared to the steady state probabihty of being at 
the separatrix n^, k{n^ — > n ~ peak) = P{nl^)/pi,{n\). 
p%{x) = X)n=^ P{n). The escape rate turns out to be: 

(8) 

We can approximate ^ in the strictly adiabatic hmit 
as for the kon calculation. Then the strictly adiabatic 
escape rate becomes: 

L [n — n ) 

koff{n^)^h{K-'')e ^<'— " " (9) 

where Umax is the number of proteins at the max- 
imum of which scales as rimax ^= V\/K^. 
The pre-exponential fatcor has the form f^iK'^'') k, 

kV/{2n)inT^ - {K'^'iVY -2K'=iV^n^^)Vl<^/n^\ More 
explicitly the escape rate from the on state scales as 
_ ^-a2n^-''y/(Cni^-K'"iv^f ^ ^j^g^g « 2 and ( ~ 
1/4 + ao/2 are constant numerical factors. The escape 
rate from the on state attractor exponentially increases 
with the increase of the equilibrium constant. A simple 
result is also obtained by replacing the effective produc- 
tion rate by the value of the effective production rate in 
the on state attractor g*^^^ in) « g'^^^ {n\): 



koff{n\) = /2(J^^')e"'-^»=^^<"^' (10) 

The equilibrium constant for the dressed genetic switch 
state in the strongly adiabatic limit is K'^^ = koff /kon ~ 

which sharply depends on the proteomic atmosphere. 

fr{K'"i) = \Jao{n^^ - {K'^W^f - 2K^iv}'^) j n^^ 
Pi « 2, /32 « 1/4 are numerical factors. 
In the weakly adiabatic regime the expo- 
nential term in the off escape rate becomes 
exp{-h/{2kW^/{K'">V^)(l/6{n\f - {n\f - 
g''ff{n\)/{5k){{n\f - {n\f))). So in the weakly 
adiabatic limit the equilibrium coefficient for the dressed 
genetic switch states K'~^^ — koff /kon scales as K'^^ ~ 

where the coefficients are determined by the positions of 
the on and off state attractors and are of the order of 
^1 « 0.01 and ^2 ~ 100. 

Whether the switch is nonadiabatic or adiabatic can 
be determined by comparing the mean free path to the 
size of the transition region. If Irsr/lmfp > 1 many 
crossings are required and the transition is adiabatic. If 
iTSr/hufp < 1 the system commits to the new attractor 
once it reaches the separatrix, hence the transition 
is nonadiabatic. In the strictly adiabatic regime the 
diffusion of the system is governed by protein diffusion 
induced by the birth-death process, as opposed to the 



weakly adiabatic regime, where diffusion due to churns 
dominates. A phase diagram showing the different 
escape mechanisms in parameter space for fixed K^'^V^ 
is shown in Figure |31 

Comparison with Numerically Exact Results 

While the mechanism of spontaneous switching or 
epigenetic escape is different in the various regimes, we 
understand the rates in all regimes using the notion of 
a stochastic separatrix. We can compare these approxi- 
mations with numerical calculations due to Kepler and 
ElstonS2i^ and our own full numerical results. 
In the nonadiabatic limit (small k — hg'^/{2k^)) the 
escape process is determined by the rate of DNA state 
fluctuations. In this regime the rates are given by 
equations ^ and |3 (gray dashed line) (Figure . These 
agree with the discrete numerical calculation of the 
mean free passage time from each basin. Our numerical 
calculations confirm that only in the extremely adiabatic 
limit (large k - flat escape rate) can the DNA fluctuations 
safely be neglected. Only for this extreme limit does 
the lifetime become determined by protein synthesis/ 
degradation fluctuations alone (light gray dashed line). 
Estimates of the input parameter would suggest that 
the weakly adiabatic regime is common for biological 
switches. In the weakly adiabatic regime the escape rate 
does not just depend on occupancy averaged growth 
rates, but still depends on the adiabaticity parameter, 
as shown in Figure ^ Neglecting the influence of DNA 
fluctuations in this limit, as many treatments have done 
would give the extreme adiabatic asymptotic value of 
the escape rate also pictured on the graph. Both the 
strictly and weakly adiabatic regimes can be obtained 
from the more general calculation using the full diffusion 
coefficient. The full treatment is only required in a small 
crossover regime (gray solid line). 
Summary 

Spontaneous transitions between attractors of genetic 
systems are caused by coupled stochastic fluctuations in 
the DNA state and protein number. Even in parameter 
regimes where the DNA state locally would appear to 
reach a steady state much more rapidly than the protein 
number state, the fluctuations due to binding and 
unbinding of transcription factors greatly influence the 
protein number fluctuations and hence modify the rate 
of spontaneous transitions between epigenetic states. 
We call such a regime the weakly adiabatic by contrast 
to the strongly adiabatic limit, where the DNA binding 
state may be taken to be in equilibrium. The mechanism 
of spontaneous switching between stable attractors in 
the weakly adiabatic regime is graphically explained by 
a churning process, which causes protein numbers to 
fluctuate from the mean flow. How the escape rates 
kon and koff depend on molecular parameters in the 
nonadiabatic, weakly and strongly adiabatic should 
allow one to understand the evolutionary constraints 
necessary to achieve stable yet responsive switches, a 
topic we hope to return to. By considering both the 
DNA and protein degrees of freedom, the rate theories 
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we have presented provide an intuitive description of 
spontaneous switching events, in terms of the molecular 
parameters that determine the functioning of a genetic 
switch. 
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